M estimators in R

manual_M_location <-
  function(x,
           k,
           type = c("Huber", "Tukey"),
           tol = sqrt(.Machine$double.eps),
           itermax = 1000) {
    
    
    mu <- mu_old <- mu_vec <- median(x) # initial value
    rob_scale <- mad(x) # it will be kept fixed
    crit <- TRUE
    iter <- 0
    
    weigth_f <- switch (type,
                        "Huber" = function(x) pmin(1, k/abs(x)),
                        "Tukey" = function(x) ((1-(x/k)^2)^2)*(abs(x)<=k)
    )
    
    while(crit){
      w_i <- weigth_f((x-mu)/rob_scale)
      mu <- weighted.mean(x = x,w = w_i)
      err <- abs(mu-mu_old)
      
      mu_old <- mu
      mu_vec <- c(mu_vec,mu)
      iter <- iter+1
      
      crit <- (err > tol & iter < itermax)
    }
    list(mu=mu, s=rob_scale, w_i=w_i, it=iter)
  }